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The effects of the rotation stress mechanism and the penetration of short- 
wave radiation below the sea-surface are examined in determining a one- 
dimensional equilibrium mixed layer depth. Starting with the Obukhov-scale 
equilibrium theory for the surface ocean boundary layer, a revised 
equilibrium theory, which includes rotation stress and radiation effects, is 
presented. This new theory is applied using climatological boundary 
conditions for the tropical Atlantic, and the results are compared with the 
observed climatological mixed layer depth. 

In general , the response of the revised model is an improvement over the 
Obukhov theory alone. Because the quality of the results are limited by 
uncertainties in the boundary conditions, no detailed evaluation of the model 
response is justifiable. However, it is concluded here that the physical 
mechanisms of rotation stress and penetration of radiation are important in 
determining a steady-state equilibrium depth of mixing for the tropical 
Atlantic. 
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I. INTRODUCTION 



The study of the oceanic turbulent boundary layer is a relatively recent 
field in Physical Oceanography. Since the pioneering work of Kraus and 
Turner (1967), a significant effort has been made in understanding and 
modelling the dynamics of turbulent mixing in the oceans. However, the basic 
contributions made by those authors still provide the theoretical foundation 
for much of the current research. 

The oceanic boundary layer or mixed layer is considered in this study to 
be the fully turbulent region that is bounded above by the air-ocean 
interface, and it is where the temperature and salinity are usually observed 
to be fairly well mixed. Below, the mixed layer is assumed to be bounded by a 
dynamically stable thermocline. 

The study of the top few tens of meters of the ocean is of considerable 
scientific interest. Ekman pumping effects, originating in the mixed layer by 
the action of wind forcing, considerably influences the dynamics of the lower 
levels and the interior circulation in general. Also, the interaction between 
the oceanic and atmospheric boundary layers is an essential mechanism 
which must be considered when making medium and long-range weather 
forecasts, since a large part of the atmospheric energy supply comes from 
the heat exchanged between these layers . Additionally, this is the region of 
primary biological productivity, which is of significant economic importance. 
An important military application is, on the other hand, the modelling of 
accoustic propagation in the oceans. 



8 



Most of the physics behind existing models of the oceanic mixed layer are 
based on the flux form of the Navler-Stokes equations of motion, with the 
Bousslnesq approximation. One-dimensional versions, like the one which is 
treated In this study, further assume the horizontal gradients of the mean 
fields to be negligible. For a steady-state situation, this will lead to the 
basic energy balance, as stated by Kraus and Turner (1967), between the 
work done by the wind stress and the surface heating, in defining an 
equilibrium depth of mixing. However, a usually neglected term in the 
turbulent kinetic energy budgets, the planetary rotation term, was recently 
examined by Garwood eta/. (1985a). They suggest a new formulation for the 
equilibrium mixed layer depth, which basically describes the Interaction 
between the zonal Reynolds stresses and the northward component of the 
planetary rotation. This new formulation might partially explain the zonal 
variation of the observed mixed layer depth In the Tropical Pacific and 
Atlantic Oceans. 

The goal of this study is to apply the Kraus and Turner (1967) and the 
Garwood et af. (1985b) one-dimensional, steady-state models on a relatively 
broad grid with realistic boundary conditions in the Equatorial and Tropical 
Atlantic Ocean, and to compare the resulting diagnostic mixed layer depths 
with the observed mixed layer depths. In this way, Insight may be gained Into 
the relative Importance of the physical mechanisms Involved. Also, we 
should be able to identify the regions and time of the year where the turbulent 
boundary layer can be represented realistically by such simple steady-state 
equilibrium models. 
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II. THE DYNAMICS OF MIXING 



A. THE TURBULENT KINETIC ENERGY BUDGET 

In the presence of a constant wind stress and a downward surface heat 
flux, the steady-state, one-dimentional turbulent kinetic energy (TKE) 
budget is, assuming no vertical advection : 

0=d/dt [E/2]= 

—u'w' du/dz -v'w' dv/dz + b'w‘-d/dz [w'(E/2)+ w'p'/p] - d, (2. 1 ) 

where E/2 is the TKE, -u'w' du/dz -v'w' dv/dz represents the production of 
kinetic energy associated with the vertical shear in the horizontal flow 
induced by the wind stress, b'w' is the buoyancy flux, -d/dz [w'(E/2)+w'p'/p] 
is the turbulent diffusion and d is the dissipation. With no loss of generality, 
we shall assume the x-axis to be oriented downstream, so that v=0. 

1 . Shear Production of Turbulence 

Let ii(z) be some monotonic function of z between the surface , z=0, 
and the bottom of the mixed layer, z=-h. For the sake of simplicity, it is 
further assumed that u(-h)=0, since no entrainement is allowed. 

Consider now (Fig. la) that some particle initially at 1 is displaced to 
some position 2 by means of a positive turbulent vertical velocity w'. 
Because du/dzz.0, the particle will acquire a negative turbulent horizontal 
velocity u' and, then, u'w' will be less tha* zero at 2. Alternatively, a 
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particle Initially at 1 , which Is displaced to 3 by means of a negative vertical 
turbulent velocity w', will acquire a positive horizontal turbulent velocity u‘ 
and, again, u’w' will be negative. In both situations, the shear production 
term -u'w' du/dz will be positive, and when this physical process is 
averaged over the whole turbulent layer, it will tend to increase the TKE at 
the expense of the mean kinetic energy. Formally, 

-uV du/dz > 0. 

Evaluated over the mixed layer, the shear production will be: 

o 

G = } -u'w‘ du/dz dz = <-u'w*> [u(0)-u(-h)], (2.Z) 

-h 

where < > represents a vertical average. At the bottom of the mixed layer, 
and under the assumption of no entrainement, -iTw'(-h) is zero. Then, 
disregarding nondlmenslonal constants, <-i7w'> Is of order -i7w'(0). Also, 
we have already assumed u(-h) to be zero. Thus, 

G = -uV(0) u(0) . (2.3) 

Furthermore, the Reynold stress -iFw’(O) can be written 



| -u’w' (0)| = |t|/p , 



(2.4) 



where p is a typical sea water density and t is the surface wind stress. 
Finally, (2.3) is rewritten as 



G = u(0)|i|/p . 



(2.5) 



2. The Buoyancy Flux 

.Under the Influence of surface heating, there will be a tendency for 
stratification, and mixing will only occur in the presence of wind stirring. In 
this case, the less dense particles near the surface will mix with the 








Fig. 1. (a) Shear production of turbulent kinetic energy, In the presence of 
wind stress, (b) Buoyant damping of turbulent kinetic energy, in the 
presence of surface heating. 
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underlying water only at the expense of the wind-generated turbulence, thus 
converting TKE Into potential energy (Fig. lb). This means that, on the 
average, the buoyancy flux -bV will be a sink term in the TKE budget. The 
equilibrium between the effects of wind mixing and surface heating Is Indeed 
the fundamental energy balance for any steady-state, one-dlmentional mixed 

layer model, as long as advectlon Is neglected. 

To relate the buoyancy flux to the surface heat flux Q 0 we apply the 

linearized equation of state and obtain: 

b'w* = agT'w', (2.6) 

where a is the thermal expansion coefficient, g Is the acceleration of 
gravity and TV Is the temperature flux. It can further be shown that, 
assuming vertical homogeneity In the mean temperature T, T'w‘ will be linear 
over the mixed layer, and 



TV(0) = -Q 0 /( pC p ), (2.7) 

where Q 0 , the surface heat flux is positive downward. 

Integrating the buoyancy flux over the mixed layer, 

© o 

B = J bV dz = ag ( TV dz = agh [TV(0)+TV(-h)]/2, (2.8) 

-h -h 

where h is the mixed layer depth. Because no entrainment Is assumed, 
TV(-h)=0. Thus, 
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B - cigh Q 0 /(2pCp) . 



(2.9) 



This simple model for the vertically Integrated buoyancy flux also 
requires the assumption that the net solar radiation be completly absorbed at 
the surface. We will see in the following sections that this assumption is 
incorrect when the vertical scale of penetration of radiation below the 
surface is comparable to the vertical scale of turbulent mixing. The effect of 
the penetration of radiation is expected to be relevant in optically clear 
oceanic waters and in summertime conditions, when the mixed layer is 
relatively shallow. In those situations, the net downward buoyancy flux may 
be significantly reduced, causing the resultant depth of mixing to be 
significantly increased. 

3. Turbulent Diffusion 

Although the diffusion of turbulence is an essential mechanism in 
transporting TKE downwards from the upper levels, it does not constitute, by 
itself, a source or sink of turbulence, which means that its net value over 
the mixed layer must vanish. Formally, 

J -d/dz lw'(E/2)+w’p7p] dz = 0. (2.10) 

-H 



Diss i p ation 

Early mixed layer models, including the Kraus and Turner (1967) 
model, recognized the existence of dissipation. However, they considered its 
influence to be negligible, given the relatively large length scale of turbulent 
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mixing. Indeed, a less careful scale analysis of the TKE budget seems to 
support the approximation. More recent models, however, include 
dissipation In the mechanical energy budget, either as a fraction of the wind- 
-stress production or as a fraction of the net TKE production. Although 
viscous dissipation may be neglected when dealing with turbulent eddies of 
length scales of order of the mixed layer depth, the dynamic interaction 
between those eddies produces other eddies of smaller and smaller 
dimensions, where the role of dissipation may finally dominate. This fact is 
easily understood with the analysis of a spectrum of turbulence, where it can 
be verified that the mechanical energy drops rapidly above a limiting 
frequency (Tennekes and Lumley, 1962). 

If such a dissipation mechanism is included in the TKE budget by 
parameterization, for example as a fraction of the TKE production, we will 
see in the following sections that the Kraus and Turner model is still usable, 
since no aditional degree of freedom is required. That was the approach of 
Geisler and Kraus ( 1969) , Niiler (1975) and others. 

For now, we will just define the vertically integrated dissipation D as: 

D = J d dz . (2.11) 

-h 



B. EQUILIBRIUM DEPTH OF MIXING 

From the results of the previous section, the vertically integrated TKE 
budget can be written as 
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B + G - D = 0, 



( 2 . 12 ) 



where B = - agh Q 0 /(p C p ) 

and 6 = u(0)|x|/p . 
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Fig. 2. The one-dimensional mixed layer equilibrium. The ‘redistribution* 
effects include the pressure redistribution and the rotation stress 
mechanism. 



16 



Equation (2.12) illustrates the basic steady-state balance between the work 
done by the wind stress and the buoyant damping (Fig. 2). Note that, because 
of turbulent mixing, the mean flow will tend to be vertically homogeneous, 
except near the surface, where turbulence Is generated, and near the bottom 
of the mixed layer, where It falls to zero. Also, the mean temperature 
profile will tend to be homogeneous. Thus, the mixed layer potential energy 
will approach its theoretical maximum value. 

In Fig. 2 an adltlonal mechanism Is Illustrated. This processs, called 
‘redistribution’, does not appear in the TKE budget, since its constituent 
terms sum to zero when the three spatial components of the kinetic energy 
budget are added. Redistribution consist of two parts: the pressure 
redistribution and the planetary rotation. Their roles are to exchange TKE 
among the x-y-z directions, tending to restore the system's Isotropy. In 
section III, the planetary rotation term will be shown to have a fundamental 
role In determining the equilibrium mixed layer depth for the Garwood et al. 
(1985a) steady-state model. 

1. The Obukhov Length Scale 

We have seen that, if dissipation is parameterized in terms of the 
other basic turbulent variables, the steadiness of the mixed layer depth will 
depend on the surface boundary conditions, i. e. , the surface wind stress and 

the surface heat flux. This Is equivalent to setting G/B = constant. 

We will now define a velocity scale u* , which is representative of the 

typical turbulent speed of the largest eddies in a wind-driven mixed layer: 

u* 2 = |i|/p . (2.13) 
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Because u(0) depends linearly on the wind speed and, thus, on the surface 
wind stress, the shear production term G in equation (2.5) can be scaled as 

G~u* 3 . (2.14) 

Combining (2.9) and (2.12-14), and setting -G/B =1, the Obukhov 
length scale is defined as 



L = h = 2 p C p u* 3 / (a g Q 0 ) . (2.15) 

Physically, L represents the scale of maximum depth of turbulent mixing, 
given the surface boundary conditions of wind stress and heat flux, for an 
equilibrium mixed layer with no advection present. 

2. The Kraus and Turner Steady-State Model 

Based on equation (2.12), Kraus and Turner (1967) considered the 
dissipation to be negligible, so that 



B + G = 0. 



(2.16) 



Also, the vertically Integrated shear production G can be written 

G= u* 3 = (|t|/p) 3/2 (2.17) 

or G= (p a C D l|2/p) 3 /2 , (2.18) 
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where Cp is the drag coefficient, p a is the air density and U is the wind 

speed measured at 10 meters above the surface. Substituting (2.9) and 
(2.18) into (2.16), the equilibrium depth of mixing is written 

h = C] ( p C p / a g Q 0 )( P a C D u2/p)3/2 , (2.19) 

where C| is a dimensionless constant of order 1 . 

More recent versions of the Kraus and Turner model assume the 
dissipation D to be proportional to the total TKE production B+G, so that 
equation (2.12) is written as 



(K + 1 )( B + G ) = 0, (2.20) 

where K is a dimensionless constant. As stated earlier, this will lead to the 
same basic result of equation (2. 19), since no aditional degree of freedom is 
included. 
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III. THE REVISED STEADY-STATE MODEL 



As noted in section I, the dynamics of the oceanic mixed layer in the 
equatorial and tropical regions are not fully understood. One of the most 
apparent but least understood features is the zonal variation of the mixed 
layer depth. 

Some of the discrepancies between the observed mixed layer depths and 
the predictions made by existing steady-state and one-dimensional models 
might be explained by the exceptionally strong vertical and horizontal 
advection which are characteristic of the tropical oceans. Also, the effect of 
penetration of short wave radiation below the surface is believed to have 
some importance for relatively shallow mixed layers, especially in the 
presence of large net solar irradiance. Another physical mechanism, 
examined by Garwood et a/. (1985a), involves the interaction between the 
zonal surface wind stress and the planetary rotation, and it may partially 
explain the zonal dependence of the steady-state mixed layer depth. 

It is beyond the scope of this work to formally treat the effects of 
advection on the dynamics of the mixed layer, and we shall restrict our study 
to steady-state situations. However, the other two mechanisms referred to 
(penetration of radiation and planetary rotation) will be included in the model 
derived in the following sections. 
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A. EFFECT OF THE PENETRATING RADIATION ON THE BUOYANCY FLUX 

The simplest model for the penetration of radiation below the sea-surface 
assumes that the net long-wave solar irradiance ( Infra-red and red ) is 
absorbed at the surface, and the short-wave irradiance decays exponentially 
with depth. Recent models ( Paulson and Simpson, 1977 ; Simpson and Dickey, 
1980) also allow for the penetration of long-wave radiation. In both 
approaches, the predicted vertical profile of absorbed Irradiance is sensitive 
to the values of two empirical parameters: the vertical scale of penetration, 
which mainly depends on the turbidity of the water, and the short-wave 
(blue-green) fraction of the net solar irradiance incident on the surface. 
Because there is only a limited amount of data on the optical properties of the 
oceans, the values of these empirical parameters will either have to be 
assumed or adjusted as tunning constants. For that reason, and also because 
a steady-state model will be applied to large oceanic regions having variable 
optical properties, the simplest model will be used, i . , e . , only the short- 
wave radiation will be assumed to penetrate below the surface. The vertical 
profile of the penetrating irradiance will have therefore the following 
exponential form 



Q 2 (z) = RQ S e z/ ^ . (3.1) 

Here, the downward short-wave irradiance, Q 2 (z), is the radiant flux density 

on a horizontal surface due to contributions from the entire upward 
hemisphere; Q s is the net solar irradiance at the surface, R is the short-wave 

fraction of the net solar irradiance ( blue-green ), z is the vertical 
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coordinate, positive upward; and X is the penetration scale, assumed 
constant with depth. 

Alternatively, the long-wave solar irradiance, assumed to be totally 
absorbed at the surface, is given by 



Qj - Q 0 R Qg , (3.2) 

where Q 0 = Q s - Q^- G>h - Q e ’ s tr,e ne * solar irradiance at the surface, minus 

the long-wave back radiation, minus the sensible heat flux, minus the latent 
heat flux. 

Assuming horizontal homogeneity of the variables, it can be shown that 
the effective heat flux Q contributing to the temperature flux over the mixed 
layer will have the form ( Garwood, 1977 ) 

Q = Q|- f l dQ 2 / dz - 2/h f dQ 2 /ds ds ] dz (3.3) 

-h r 

or Q = Q 0 - R Q s [ 2\/h - e ~ h/ ^ ( 1 + 2\/h )]. (3.4) 

For deep mixed layers ( h/\ >>1) or, puting R = 0, expression (3.4) 
reduces to Q = Q 0 , which corresponds to the net heat flux being absorbed at 

the surface. For shallow mixed layers ( h/\<<l), expression (3.4) reduces 
to Q = Q 0 - R Q s . Under these circumstancies the effect of the penetrating 

radiation on the effective heat flux is maximum. 

It is also apparent from expression (3.4) that, in summertime conditions 
( small h, large Q s ) and in clear oceanic waters ( large X), Q can be 
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significantly smaller than Q 0 , thus reducing the buoyancy flux and increasing 

the vertical extent of the wind generated turbulence. 

From expression (2.9), the effective vertically integrated buoyancy flux 
will then be 



B = - ag h Q /( 2pC p ), 



(3.5) 



where Q is given by (3.4). 



A. THE ROTATION STRESS MECHANISM 

The steady-state turbulent kinetic energy budgets for the three spatial 
components ( x eastward, y northward, z upward ), vertically integrated 
over the mixed layer are: 

0 = G x ^ P x + ^ u v (0) ^ Oy ht x / p — D/3 (3.6a) 

0 = Gy + Py-Q z h uV(0) -D/3 (3.6b) 

0 = B P 2 — Oy hT x /p — D/3. (3.6c) 

Here G x and Gy are the horizontal components of the shear production, B Is 
the vertically integrated buoyancy flux, the Pfs are the components of the 
vertically Integrated turbulent transport, and t x Is the zonal component of 
the surface wind-stress. The Qi's are the two spatial components of the 

planetary rotation and D is the dissipation, assumed to be equally partitioned 
among components (local isotropy assumed). Since the transport and 
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planetary rotation terms vanish, the sum of the three above equations simply 
leads to : 



B + G - D = 0. (2.12) 

An aditional process, not apparent in (2.12), is evident in equations 

(3.6). This is the exchange between horizontal and vertical turbulent kinetic 
energy, in the presence of the zonal wind stress t x and the meridional 

component fi y of the planetary rotation. Note that the sign of t x will 
determine whether the term Qyht x /p is a source or sink of horizontal 
turbulent kinetic energy. With easterly winds ( t x < 0 ), the horizontal 

kinetic energy is expected to be enhanced at the expense of the vertical 

kinetic energy, with a consequent shallowing of the mixed layer. With 
westerly winds ( t x > 0 ), the exchange between the vertical and horizontal 

components will have the opposite sign, leading to an increase of the vertical 
extent of turbulent mixing. This mechanism, usually overlooked in planetary 
boundary layer models, was first examined by Garwood et al. (1985b), in an 
attempt to explain the unusually deep extent of mixing in the central and 
western Pacific. 

If the dissipation is parameterized in terms of the vertical average of the 
total turbulent kinetic energy, i.e, 

0=m!<E> 3/2 , (3.7) 
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where mi and subsequent m^'s are dimensionless constants, ( 2 . 12 ) can be 
written 



G + B - mi < E > 3/2 = 0 . (3.8) 

Also, the vertical component P 2 of the turbulent transport can be 

parameterized in terms of the vertical average of the TKE. According to Rotta 
( 1951 ) , a first order approximation to the transport term P z is 



P 2 = m 2 <E> 3/2 . (3.9) 

Substituting (3.9) into (3.6c) and combining with (3.8), will yield a 
system of two equations which can be solved for the mixed layer depth h and 
the vertical average of the turbulent kinetic energy < E >, given the surface 
boundary conditions and the meridional component of the planetary rotation. 
The predicted value of h will then reflect the interaction between the zonal 
wind stress and the Earth rotation. 

C. EQUILIBRIUM DEPTH OF MIXING 
We define 





L = 2pC p u* 3 /( ag Q ) 


(3.10) 




G = m 3 u * 3 


(3.11) 


and 


$ = Oy^x^p/ ( QQ Q)> 


(3.12) 
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where the dimensionless parameter $ is the quotient between the conversion 
of vertical turbulent kinetic energy by means of the rotation mechanism and 
the buoyant damping of vertical turbulent kinetic energy. Then, a diagnostic 
equation for the mixed layer depth can be derived: 

h = q L /( 1 + C 2 *), (3.13) 

where the non-dimensional constants and C 2 are related to the m{s: 

= m3(m 2 /m 1 -1/3)/(m 2 /m 1 +2/3) (3.14) 

and C 2 = 2/(m 2 /m 1 +2/3). (3.15) 

Note that the effective surface heat flux Q in equations (3. 10) and (3. 12) 
reflects the penetration of radiation below the surface, as modelled by 
equation (3.4). Thus the diagnostic value of h will be dependent on the values 
of four tunning parameters: the C(s, the factor R and the vertical scale of 

penetration X. Because Q is a relatively complex function of h, a solution for 
h from (3.4) and (3.13) is only possible by iterative techniques. An 
alternative method, which significantly reduces the computational effort, is 
to substitute the Obukhov length scale L for h in equation (3.4), and then 
solve for Q. This simplification is justified and will be utilized in the next 
sections, since the primary goal of this study is to examine the gross 
characteristics of the solution space. 

Substituting (3.13) back into (3.7), the following expression for the 
vertically integrated dissipation is obtained: 
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D = m t < E > 3/2 = C 2 m 3 u* 3 ( 1+2*)/l2( 1+C 2 * )]. 



(3.16) 



Physically, the dissipation must always be non-negative, which yields 



1+2* 2 0 

or *2-1/2. (3.17) 




Fig. 3. The dimensionless depth scale h* plotted against the rotation stress 
variable <D . 



27 



For <t> < -1/2, and in absence of advection, no steady-state situation is 
possible and the mixed layer will continually entrain. Also, noting that, in 

equation (3.13), the denominator must be positive, the limiting case 
expressed by (3.17) suggests that C£ < 2 for all situations (since C;? is 

assumed to be a constant). Zeman and Tennekes (1975) obtained 
m2/mi~1/2, which gives 0% ~ 12/7. Using the suggested value, a 

nondlmenslonal mixed layer depth Is given by ( Fig. 3 ) 

h* = h/ C^L = 1/( 1 + 12*/ 7 ). (3.18) 

Note that to the limiting value ♦ = -1/2 corresponds a theoretical maximum of 
h* = 7. 
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IV. THE EQUILIBRIUM MIXED LAYER DEPTH IN THE TROPICAL ATLANTIC 



A. THE BOUNDARY CONDITIONS 

In the following sections, the method of application of the steady-state 
mixed layer model to an oceanic macro-scale domain is presented. The goal 
of this study is to gain insight into the relative importance of the physical 
mechanisms involved, as well as to identify the regions and the time of the 
year where the observed turbulent layer can be represented realistically by 
such a simple steady-state model. This goal is limited to some degree by the 
availability and suitability of experimental observations of the boundary 
conditions. Although a significant effort has been made in the last decade to 
obtain accurate and detailed measurements of the surface boundary 
conditions, such data are frequently confined to limited oceanic areas. The 
available climatological data for large oceanic areas are often unreliable, 
sparse and distributed irregularly, depending on the usual routes followed by 
mariners. Given these limitations, the tropical Atlantic Ocean is the best 
covered ocean basin. 

Except for the observed surface wind-stress, where the digitized data of 
Hellerman and Rosenstein (1983) were used, all the other data were 
interpolated from climatological atlases. The area chosen was the tropical 
and equatorial Atlantic between 6 degrees south and 30 degrees north. 

1 . The Observed Mixed Laver Depth 

The atlas of Robinson et at. (1979) was used to obtain the observed 
monthly values of the mean depths to the top of the thermocline for the 
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equatorial and tropical Atlantic. The top of the thermocline is therein defined 
as the depth at which the temperature is 1 . 1 degrees centigrade less than the 
surface temperature. Although this definition differs from the usual 
definition of 'mixed layer depth' (the vertical extent of turbulent mixing), the 
difference should not be significant for the present study, provided density is 
mostly dependent upon temperature, and the boundary layer is well mixed. 
For the cases where the mixed layer mean temperature decreases 
significantly with depth, the above criterion might underestimate the vertical 
extent of mixing. 




Fig. 4. Zonal variation of the mean mixed layer depth for the month of 
September (equator). The curve is the result of an adjusted cubic spline fit 
over the observed values, represented in the figure by circles. 
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Because the contours displayed in the atlas have a relatively large 
increment of 15 meters, a cubic spline fit was made for each parallel of 
latitude, and the resulting curve was further adjusted so that the interpolated 
values for each interval would lie between extremes. The resulting maximum 
graphical error is estimated to be ±10 meters. No indication is given by the 
authors about the estimated errors for their analized depths. 

Fig. 4 shows the zonal equatorial profile of the mean mixed layer 
depth for the month of September, as computed by the above technique. 

2. The Surface Heat Flux 

The heat budget atlas of Hastenrath and Lamb (1978) was used to 
provide best estimates of the net surface heat flux and the net surface solar 
irradiance for the area. For the assumptions made and the calculation 
details, the reader is referred to the introductory part of the atlas. 

The contour interval for the surface heat flux is 40 W/m 2 , which 
constitutes a low resolution. In order to minimize the resulting errors, a 
cubic spline surface was first fitted over a grid of digitized values, and the 
resulting interpolated values were further adjusted for each parallel of 
latitude. Thus, for each interpolation interval, the calculated values would 
lie between extremes. The considerable effort required for this process is 

well Justified because of the sensitivity of the diagnostic equations to small 
deviations in the values of Q 0 , specially when Q 0 is small. 

The relative and absolute errors in the net surface heat flux, as 
indicated by the authors, are estimated at less than 10 and more than 
20 W/m 2 , respectively. The absolute error in interpolated values is 
estimated to be less than 20 W/m 2 . 
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As an example, Fig. 7a shows the zonal equatorial profile of the mean 
surface heat flux for the month of September, as determined by the above 
technique. 

3. The Surface Wind Stress 

The surface wind stress data are those of Hellerman and Rosenstein 
(1983), for a mesh grid of 2 degrees of latitude by 2 degrees of longitude, 
which is suitable for our purpose. The two components are calculated by 
means of the bulk aerodynamic formulation: 

T x = PaCoU x (U x 2 + U y 2 )1/2 (4.1) 

and T y = p a Co U y (U x 2 + Uy 2 )l/2 # (4.2) 

where U x and Uy are the two components of the wind speed at 10 meters 
above the sea-surface. The drag coefficient, Cp, depends upon wind speed 

and stability, according to the formulation of Bunker (1976). The maximum 
standard error in the surface wind-stress values for the area considered is 
estimated to be less than 0.01 N/m 2 . 

The wind stress data are given for each odd whole degree of latitude 
and longitude, and in the present scheme we will be working with even 
parallels of latitude. Thus the original data had to be interpolated. The 
resulting interpolation errors are considered to be negligible for the purpose 
of the present study. 

4. Steady-State versus Time-Dependent Solutions 

Formally, a steady-state mixed layer model should only be applied to 
a real oceanic situation when entrainment has stopped and the top of the 
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thermocllne has retreated to an equilibrium level, dependent on the balance 
between the downward surface heat flux and the surface wind-stress. On a 
seasonal time scale, l.e., neglecting daily variations, and for mid- 
latitudes, the above situation is expected to occur In spring and summer, 
when the mean downward surface heat flux is rising and the wind speed is 
decreasing. For the equatorial and tropical regions, however, this process 
is complicated, since the mean surface heat flux follows a complex annual 
variation, while the wind speed has a strong annual component. 

For simplicity, the diagnostic mixed layer depth should be calculated 
for the same month throughout the chosen area, though the requirement of 
steadiness will not be satisfied In some regions. A logical and simple way to 
solve this problem is to plot the annual time series of the zonal profiles of 
the observed mixed layer depths for some parallels of latitude and choose 
the time of the year which best conforms to the condition of steadiness for all 
latitudes. Figure 5 shows these plots, as compiled and further interpolated 
from Robinsons's atlas. The month of September was chosen as conforming 
best to these conditions. Except at the equator, where the top of the 
thermocline seems to be still deepening, the other latitudes show 
approximate steady-state situations, during September. 

The analysis of Fig. 5a reveals a strong annual component of the 
equatorial mixed layer depth, which suggests the domination of the wind 
stress In determining the vertical extent of mixing. In Fig. 6, the annual time 
series of the net surface heat flux and the wind speed are shown. Indeed, a 
qualitative compairason of the contour lines of both figures (5a and 6) 
reveals a much stronger correlation between the mixed layer depth and the 



33 




■ 


> 90 m 


mi 


75-90 


l/.v-HIT 

12222 


60-75 


□ 


45-60 


ED 


30-45 


□ 


15-30 


n 


< 15 



Fig. 5. Annual zonal variation of the mean mixed layer depth, (a) Equator. 
( b) 5 degrees north. 
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Ftg. 5(cont. ). Annual zona) variation of the mean mixed layer depth . (c) 10 
degrees north, (d) 15 degrees north. 
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Fig. 6. (a) Annual zonal variation, in W/m 2 , of the net surface heat flux at 

the equator, (b) Annual zonal variation, in meters, of the wind speed at the 
equator. 
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wind speed than between the mixed layer depth and the surface heat flux. In 
particular, the wind speed maximum (August to October) and minimum 
(March to April), coincide well with the corresponding relative extremes in 
the mixed layer depth. Similar agreement is not apparent between the mixed 
layer depth and the surface heat flux, as so far as an equilibrium state is 
concerned. 

On the other hand, one might expect the annual variation of the surface 
heat flux to show a marked response due to the annual variation of the Sun's 
declination, i.e., to show a strong biannual component. However, the 
exchange of heat between the ocean and the atmosphere is, in the equatorial 
and tropical regions, significantly Influenced by other variables. These 
variables, such as cloudiness and precipitation, have seasonal variations 
which are more complex and difficult to explain. 

During the months of May to September, another Interesting feature in 
the eastern part of the Equator is the mixed layer depth minimum. This 
occurs when the surface heat flux is negative (upward). Without advection, 
an equilibrium situation is obviously precluded. Interestingly, the atlas of 
Robinson et a/. (1979) shows that this is also a zone of strong equatorial 
upwelllng and surface divergence. Thus, an advective mechanism may 
largely influence the dynamics of turbulent mixing for that region. 

From the above discussion and also from the knowledge that the 
equatorial and tropical regions are zones of important oceanic heat 
divergence, it might be concluded that horizontal and vertical advection play 
an important role in determining the vertical extent of mixing at the equator. 
As lateral transport in the ocean diverts much of the heat exchanged with the 
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atmosphere, this effect would tend to spaclally homogenize the surface heat 
flux available for the vertical turbulent processes. Thus, the seasonal 
variations of the mean mixed layer depth would primarily reflect the 
seasonal variations of the wind speed. 

B. METHOD OF MODEL TUNING 

The models described in sections II and III shall be applied to a set of 
some parallels of latitude in the tropical and equatorial Atlantic, from 6 
degrees south to 30 degrees north. For each parallel, the spacing of 
calculated values will be I degree of longitude. The purpose of the following 
procedure is to identify the zones where the observed mixed layer is well 
represented by the above steady-state models, as well as to analyse the 
relative influence of the physical mechanisms discussed in section III: the 
penetration of radiation below the surface and the rotation stress 
mechanism. For each parallel of latitude, a zonal profile of the following 
nondimensional depth will be first calculated: 

H* = H (1+ C 2 0)/L , (4.3) 

where H is the observed mixed layer depth and C 2 = 12/7. To calculate the 
Obukhov length scale L and the rotation stress variable 0, the penetration 
parameters will be set to , respectively, R=0.2 and \ = 6 meters. Note 
that, for each longitude point, H* is numerically equal to the constant C^. In 

the ideal case of a perfect agreement between the model's results and the 
observations, the resulting curve of H* for each latitude would be a constant 
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of order 1. The departure of the calculated values from that Ideal solution 
will constitute a method of identifying the zones where steadiness is not 
verified or advectton plays an important role. For the situations where the 

surface heat flux is upward, 1. e. , where a one-dimensional equilibrium is not 
possible, H* will be set to 0. The described procedure also makes possible 

the selection of the areas where the models will be applied. For these areas, 
the below defined diagnostic mixed layer depths will then be calculated: 

(i) h 0 : the equilibrium mixed layer depth, assuming R=0, C2=0. This 

corresponds to the steady-state Kraus and Turner diagnostic depth, where 

radiation is assumed to be totally absorbed at the surface. 

(ii) hj : the equilibrium mixed layer depth, assuming R=0.2, X = 6m, 

C2=0. This is the same as h 0 , except that radiation is assumed to penetrate 

below the surface, according to the model given by equation (3.4). 

(iii) f»2 •' the equilibrium mixed layer depth, assuming R=0.2, X = 6m and 

C2*0. Here, the rotation stress mechanism is also included, according to the 
formulation of Garwood eta/. (1985a). 

The pertinent equations are given in sections II and III, and the following 
constants are assumed: 



p = 1025 Kg/m 3 
p a = 1.2 Kg/m 3 

a = 2x 10 -4 C _1 
Cp = 3890 JKg-1 C" 1 
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g = 9.8 m/s 2 
0 = 7.29 x 10~5 S "l. 



For each of the above defined depths and for each latitude, the tunning 
constants C| and C 2 will be adjusted, so that the summation of the squares 

of the differences between the diagnostic mixed layer depths and the observed 
mixed layer depths is a minimum, and C|i 0 and 0 i C 2 i Z. Following this 

tuning process, a standard deviation will be then computed, according to the 
definition: 



N 

s2 = 2CH i -h i ) 2 /(N-l), (4.4) 

itl 

where the Hfs are the observed mixed layer depths, the hfs are the 

diagnostic mixed layer depths and N is the number of longitude points. The 
standard deviations will serve as a quantitative verification of the models, 
the model assumptions, and the accuracy of the boundary conditions. 

As discussed earlier, the rotation-stress variable 0 is required 
theoretically to be greater than -1/2, the value corresponding to a situation 
of zero dissipation. For that reason, G) will be set equal to -1/2 whenever 
its calculated value is less than -1/2. During the month of September, this is 
likely to occur in the western Atlantic, where the zonal surface wind stress 
is large and the surface heat flux is small. 

The selection of the factor R and the penetration depth X is an 
imprecise process because of the paucity of data for the optical properties of 
the equatorial and tropical oceans. Although the choice here of values for R 
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and X. was based on a preliminary tunning of the model at a specific site, it is 
still somewhat arbitrary. Without a more complete data set, the inclusion of 
these parameters as spacially-varying tuning parameters does not seem 
justifiable at this time. Such a procedure would tend to mask the effects of 
other important physical mechanisms, such as horizontal and vertical 
advection. An alternative might be to define two different scales of 
penetration for the long and short-wave radiation, and then to use the 
experimental results of Paulson and Simpson (1977). However, this would 
also complicate the problem. Furthermore, Paulson and Simpson's data are 
valid only for a specific spatial and temporal situation in the North Pacific. In 
the absence of detailed measurements of the ocean’s optical properties, 
there is no simple solution to this problem, and the results of the present 
work will be somewhat affected. 
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C. RESULTS 

Appendix 1 shows the zonal variation of the nondimensional depth H* for 
latitudes 0, 10, 20 and 30 degrees north. Only two particular cases, the 
equator and 10 degrees north, will be presented and discussed in the text. 

1 . Equator 

Fig. 7 shows the equatorial boundary conditions : the surface heat flux 
Q 0 , the net surface solar irradiance Q s , the total surface wind stress T, and 
the zonal wind stress T x . Fig. 8a shows the zonal variation of the 
nondimensional depth H*. computed using equation (4.3). It is apparent from 

this result that only a limited zone, between about 31 degrees west and 4 
degrees west seems to be in a steady-state balance. For some of the 
remaining regions, the presence of an upward total surface heat flux 
precludes the possibility of an equilibrium, causing the top of the thermocline 
to continually entrain. We have previously suggested that the turbulent 
processes near the eastern boundary are likely to be dominated by equatorial 

upwelling, which can reduce significantly the vertical extent of mixing. The 
values from Fig. 7a seem to confirme this hypothesis, since the values of H* 

are less than 1 near the coast. For the western region, other mechanisms 
might partially explain the departure of the model's results from the 
observations. The first of those is related to the strong wind-driven current 
near the coast of Brazil, the “Guiana current". This flow and the associated 
meridional variations in the surface heat flux might cause the local dynamics 
to be dominated by horizontal advection. A second possibly contributing 
factor is the influence of the Amazon River runoff. According to the data on 
the mean surface salinities in Robinson's atlas, this runoff extends to more 
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Fig. 7. Mean surface boundary conditions at the equator, during September . 
(a) Net downward surface heat flux, (b) Net downward surface solar 
irradiance. 
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Fig. 7(cont.). Mean surface boundary conditions at the equator, during 
September, (c) Total surface wind stress, (d) Zonal surface wind stress. 
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than 300 nautical miles eastward from the coast and it affects the upper few 
tens of meters of the water column. As examined by Garvine (1984), an 
advected shallow plume of buoyant water can decrease the mixed layer 
depth. Also, since the water advected from the Amazon is more turbid than 
the average sea water, the value chosen for the penetration scale A. is 
probably too large for the region, causing the model to overestimate the 
mixed layer depth. Fig. 7a shows precisely this effect : approaching the 
western boundary, H* falls markedly to about 0.4. 

For the equatorial region between 31 degrees west and 4 degrees west, 
the diagnostic values were calculated for the mixed layer depths scales h 0 , 

hi and h 2 . The results were compared with the observed mixed layer depths. 

Fig. 8 shows the zonal profiles obtained for each of the defined scales. The 
vertical bars represent the effect of varying Q 0 by ±10 W/m2, which 

indicates how sensitive the model is to moderate errors or variations in the 
heat flux data. The comparison of Figures 7a and 8 also shows that h 2 is 

significantly more sensitive to variations in Q 0 than are the other depth 
scales, especially when Q 0 is small. In determining a best fit for h 0 and hi , 
the corresponding standard deviations were s 0 = 5.4 meters and si= 6.2 
meters. These optimized values were obtained with Cp 1 .3 and Cp 1 .0, for 
h 0 andhi, respectively. 

Fig. 9 illustrates the tunning process of Cj and C 2 , in computing the 
zonal variation of the depth scale h 2 « Displayed here are isopleths, in m2, of 

the sums of the squares of the differences between the diagnostic and the 
observed depths, as a function of Ci and C 2 values. For each grid point of 

C] and C 2 values, a zonal profile of the diagnostic depth h 2 was calculated, 
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Fig. 8. Model results for the equator, (a) Zonal variation of the 
nondimensional mixed layer depth H«. (b) Zonal variation of the Kraus and 

Turner diagnostic mixed layer depth h 0 . The circles represent observed 
values, while the small squares represent interpolated values. The vertical 
bars show the effect of varying the surface heat flux Q 0 by ±10 W/m 2 . The 

shaded areas in (a) correspond to the zones where a one-dimensional, 
steady-state balance is assumed impossible. 
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Longitude 



Fig. 8(cont.). Model results for the equator, (c) Zonal variation of the 
Kraus and Turner diagnostic mixed layer depth hj. (d) Zonal variation of the 

Garwood diagnostic mixed layer depth h 2 . The circles represent observed 
values, while the small squares represent interpolated values. The vertical 
bars show the effect of varying the surface heat flux Q 0 by ±10 W/m^. 
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as well as the zonal summation of the squares of the differences between 
those depths and the observed mixed layer depths. The arrows in Fig. 9 

indicate the optimal values of the constants, corresponding to a minimum 
standard deviation = 3.3 meters. 

There are some interesting features to be noted in the results: 

(i) There is no quantitative improvement in the quality of the results 
for the Kraus and Turner model (h 0 ) when the effect of the penetrating 

radiation is included (hj). This is explained by the fact that the surface heat 
flux is relatively large in that region, and the vertical scale X of penetration 
is much less than the computed Obukhov length scale. 

( ii ) The zonal variation of the mixed layer depth near the western limit 
is better represented when the rotation stress mechanism is included (h^). 

In that region, the zonal wind stress has relatively large negative values, 
which causes d> to also have relatively large negative values. 

(iii) The model seems much more sensitive to the value of C| than to 
the value of C?. Fig. 9 shows a relatively large zone, between C2=1.2 and 
C 2 =l .6 where the standard deviation remains practically constant. This fact, 
associated with the similarity of the results between h 0 and hj , might indicate 

that the effect of the rotation stress mechanism is small, for the region taken 
as a whole. 

In general, the revised model shows a remarkably good agreement 
with the observations. However, and because the Kraus and Turner h 0 's are 

already, per se, in good agreement with the observed values, the results 
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are not conclusive concerning the Importance of the penetrating radiation and 
the rotation stress mechanism. The eastern part of the equator, where the 
observed zonal wind stress Is positive, would constitute a suitable domain to 

further test the significance of the rotation stress mechanism. There, the 
diagnostic depths h? would be larger than the Kraus and Turner depths h 0 and 

hj. However, as we have seen, other mechanisms appear to affect the 

dynamics of mixing In the region, and will need to be included In such a test. 




Fig. 9 . Variation of the sum of the squares of the differences, in .between 
the observed equatorial mixed layer depths H and the diagnostic mixed layer 
depths h2 with the tunning parameters and C2. 
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Z. 10 degrees North 



Fig. 10 shows the mean boundary conditions for 10 degress north 

during the month of September. In Fig. 11a, the zonal variation of the 
nondimensional depth H* is displayed. This figure shows that the zone which 

seems to be in a steady-state balance (between about 38 degrees west and 17 
degrees west) is considerably smaller than the corresponding zone for the 
equator. Although the surface heat flux is always positive for the entire 
domain, the values are usually smaller than at the equator. When the effect 
of the penetration of radiation is included, the result is to further reduce the 
effective heat flux at the surface, as defined by equation (3.4). This effect 
reduces buoyant damping and precludes the possibility of an equilibrium for a 
larger area. Although unsteadiness may partially account for these results, 

advection seems to be a more plausible explanation for the regions where the 
effective heat flux is positive, since H* is always less than 1. A meridional 

gradient is present in the surface heat flux because of the strong wind-driven 
current in the western part of the tropical Atlantic, and this effect might 

cause buoyant water to be advected northward, reducing the vertical extent 
of mixing. The variation of H* near the eastern boundary is more difficult to 

explain. From about 18 degrees west, where it has a value of 1 , H* continues 

to increase eastward, which indicates that the observed depth becomes 
larger than the diagnostic depth h 2 - The analysis of the temperature 

distribution in the atlas of Robinson a/. (1979) shows that this is a region 
of particularly strong coastal upwelling. However, positive vertical 
advection can only reduce the depth of the turbulent boundary layer, as 



50 





Fig. 10. Mean surface boundary conditions at 10 degrees north, during 
September, (a) Net downward surface heat flux, (b) Net downward surface 
solar irradiance. 



51 





Fig. 10(cont). Mean surface boundary conditions at 10 degrees north, during 
September, (c) Total surface wind stress, (d) Zonal surface wind stress. 
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stated by Garwood et at. (1985b). Thus, this mechanism cannot be 
responsible for the difference. 

The diagnostic depths h 0 , h^ and hi were calculated for the region 

between 38 degrees west and 17 degrees west. Fig. 11 shows their zonal 
variations. The vertical bars represent the effect of varying Qo by ±5 W/m 2 , 
from which it can be concluded that the model is much more sensitive here to 
small variations in the heat flux than it was for the equator, especially when 
the penetration of radiation effect is included. This is explained by the 
already discussed fact than the surface heat flux at 10 degrees north is 
significantly smaller than it is at the equator. In some areas, particularly in 
the western half of the domain, the sensitivity is so large that the resulting 
depths vary by a factor of three when those effects are considered. 
Considering that the relative errors in the surface heat flux data, as 
estimated by Hastenrath and Lamb (1978), are less than 10 W/m 2 , only a 
qualitative interpretation of the results for this latitude should be made. 

For each of the diagnostic scales h 0 , h^ and h 2 , the optimized tuning 

constants and the corresponding standard deviations have the following 
values: 



h 0 : Cj = 2.6, C 2 = 0.0, s = 11.8 m 
hj : Cj « 0.8, C 2 = 0.0, 3 = 7.1 m 
h 2 : Cj = 0.8, C 2 = 0.7, s = 5.9 m. 



A significant improvement in the quantitative quality of the results is 
apparent, when the penetration effect is included (hi ). Note that the optimum 
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Fig. 11. Model results for 10 degreees north, (a) Zonal variation of the 
nondimensional mixed layer depth H«. (b) Zonal variation of the Kraus and 

Turner diagnostic mixed layer depth h 0 . The circles represent observed 

values, while the small squares represent interpolated values. The vertical 
bars show the effect of varying the surface heat flux Q 0 by ±5 W/mZ. The 

shaded areas in (a) correspond to the zones where a one-dimensional, 
steady-state balance is assumed impossible. 
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Fig. 1 l(cont). Model results for 10 degrees north, (c) Zonal variation of the 
Kraus and Turner diagnostic mixed layer depth hp (d) Zonal variation of the 

Garwood diagnostic mixed layer depth h 2 . The circles represent observed 
values, while the small squares represent interpolated values. The vertical 
bars show the effect of varying the surface heat flux Q 0 by ±5 W/m^. 
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tuned value of Cl was reduced from 2.6 to a more plausible 0.8, which 
approaches the value found for the equator. In general, the inclusion of the 
penetration effect tended to homogenize the diagnostic mixed layer depth 
field. This caused a better agreement between the results and the 
observations, particularly near the extremes of the longitudinal interval. As 
expected from the analysis of the boundary conditions for the zonal wind 
stress, the western part of the domain is better represented when the 
rotation stress mechanism is included. Also, a small improvement is 
observed near 22 degrees west, where the observed mixed layer depth has a 

minimum. This is explained by the zonal wind stress being positive in that 
region. When the eastern boundary is approached, T x becomes more and 

more positive, causing the diagnosed value for h? to be less than the 
observed value. 

Fig. 11 illustrates the tuning process for h£. Again, the contour lines, 

in m 2 , of the sums of the squares of the differences between the diagnostic 
and the observed depths, are displayed as a function of and C£. The 

optimal value of Cj = 0.8 has not changed with the inclusion of the rotation 
stress mechanism. This tends to verify the merits of the procedure, in the 
sense that no other physical mechanisms, absent from the model's 
formulation, are being compensated for. On the other hand, although Fig. 11 
shows that the model is still more sensitive to variations in Cj than to 
variations in C£, the optimal value of C 2 is much better defined here than at 

the equator. Also, this relative sensitivity seems to verify the importance of 
the rotation effect. 
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Although a much better representation of the mixed layer depth field 
was obtained with the revised model at 10 degrees north, this verification Is 
limited by the uncertainty in the boundary conditions. As we have already 
noted, the model is so sensitive to small variations in the heat flux data, that 
no detailed interpretation of the results is justifiable. 




Fig. 12. Variation of the sum of the squares of the differences, in m^, between 
the observed mixed layer depths H and the diagnostic mixed layer depths h 2 

with the tunning parameters Cj and C 2 ( 10 degrees north). 
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V. CONCLUSIONS 



Starting with the theoretical frame work of the one-dimensional, 
steady-state Kraus and Turner ( 1967) model for the surface ocean boundary 
layer, a revised theory was presented. This new theory includes the effect of 
the penetrating radiation below the surface on the buoyancy flux, as well as 
the rotation stress mechanism examined by Garwood ( 1935a, b). In the 
presence of a downward surface heat flux and a surface wind stress, the 
penetration of radiation will tend to reduce the buoyant damping of 
turbulence, enhancing vertical mixing. The rotation stress mechanism, 

driven by the interaction between the meridional component of the planetary 
rotation O y and the zonal wind stress T x , predicts the equilibrium mixed 

layer depth to increase when i x is negative (westward) and to decrease when 
T x is positive (eastward). 

To test its response, and to gain insight into the relative importance of the 
physical mechanisms involved, the revised model was applied to a set of 
boundary conditions in the tropical Atlantic. Except for the surface wind 
stress, where digitized values were used, the other boundary and oceanic 
conditions (observed mixed layer depth, surface heat flux and surface 
downward solar irradiance) were interpolated from low-resolution 
climatological atlases. Since the model is particularly sensitive to small 
variations in the heat flux data, the quality of the results is considered to 
have been somewhat affected by the amount of uncertainty in the boundary 
conditions. Also, in the absence of detailed measurements of the optical 
properties of the tropical oceans, the values of the penetration parameters A. 
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and R were assumed to be constant over the whole domain. Given these 
limitations, no detailed interpretation of the results is appropriate here. 

In general, the response of the revised model compares favorably with 
the observations. The improvement of the results when the two new 
mechanisms are included is very significant, particularly at 10 degrees 
north. This was expected, in so far as the penetration effect is concerned, 
since at 10 degrees north the mixed layer is relatively shallower than at the 
equator, and the downward surface heat flux is generally smaller. For this 

latitude, an important result which seems to support the need to include the 
rotation stress formulation is that the optimal value of the constant C] 

remains constant when that process is included in the tuning of the model. 

In a recent paper, Garwood et at. (1985b) examined the zonal dependence 
of the climatological mixed layer depth in the near-equatorial Pacific and 
concluded that the rotation stress mechanism is a plausible explanation for 
the deep mixing in the central equatorial Pacific. Although that feature is not 
so pronounced in the tropical and near-equatorial Atlantic, the results 
obtained here seem to support the same conclusion, especially in its western 
part. On the other hand, the tropical Atlantic is a relatively smaller ocean 
basin, and the influence of coastal physical processes, like coastal 
upwelling, are more likely to influence significantly the dynamics of the 
boundary layer. Indeed, this effect is apparent from the results, since only 
the central part of the domain seems to be in a one-dimensional equilibrium 
state. Also, the sensitivity of the model to small variations in the surface 
heat flux shows that this equilibrium state is, as it is in the tropical Pacific, 
a delicate one. 
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The choice of the tropical Atlantic as a domain for application of the model 
was primarily motivated by the existence of a suitable set of surface 
boundary condition data. Future research on the subject should extend the 
domain to the non-tropical regions and should include unsteadiness and 
advection. Although the model derived in this paper is a simplistic one- 
-dimensional representation of a complex phenomenon, and the quality of the 
results were probably affected by the lack of detail of the available data, it 
is concluded here that the physical mechanisms of rotation stress and 
penetration of radiation are important in determining a steady-state 
equilibrium depth of turbulent mixing for the tropical Atlantic. 
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APPENDIX A - MODEL OUTPUT 





Longitude 



Fig. 13. Zonal variation of the nondimensional mixed layer depth H*. Top : 

equator. Bottom : 10 degrees north. The shaded areas correspond to the 
zones where a one-dimensional, steady-state balance is assumed 
impossible. 
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Longitude 



Fig. 14. Zonal variation of the nondimensional mixed layer depth H*. Top : 20 

degrees north. Bottom : 30 degrees north. The shaded areas correspond to 
the zones where a one-dimensional, steady-state balance is assumed 
impossible. 
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APPENDIX B - COMPUTER PROGRAMS 



' Program depth - calculates the steady-state mixed layer depth 

variation for a given parallel of latitude. 

' Joaquim Filtpe Gaspar - April 1985 - Naval Postgraduate School. 

DEFINT i-k , n 

DIM long(200) ,Qo(200) , Qs(200) , tau(200) ,taux(200) , depth (200) , h( 200) 

******************** Variables ******************************* 

' lat - latitude in degrees ; long - longitude in degrees 

' Qo - surface heat flux in W/m 2 ; Qs - surface solar irradiance in W/m 2 
' Q - effective heat flux in W/m 2 

' tau - total wind stress in N/m 2 ; taux - zonal wind stress in N/m 2 
' depth - observed mixed layer depth in meters 
‘ h - diagnostic mixed layer depth in meters 
‘ L - Obukhov length scale in meters 
‘ Fi - rotation stress parameter 
' R,lb - radiation penetration parameters 
’ cl ,c2 - model tunning parameters. 

* ************************************************************** 

Inputdata: 

INPUT ’ Enter latitude lat : 1 at$=STR$ ( 1 at ) 
f l$="sh*+lat$ : OPEN "in’, #1 ,f1$ 
f2$='qs’+lat$ : OPEN ’in‘,#2,f2$ 

f3$='tau’+lat$ : OPEN ’in’,#3,f3$ 
f4$='taux*+lat$ : OPEN 'in',*4,M$ 
f5$=’mld”+lat$ : OPEN ’in', #5,f5$ 



n=0 

WHILE NOT EOF(l) 

n=n+l 

INPUT# 1 , long(n) ,Qo(n) 
INPUT#2,lg,Qs(n) 
INPUT#3, lg,tau(n) 
INPUT#4,lg,taux(n) 
INPUT#5,lg,depth(n) 
WEND 
CLOSE 
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Inputparameters: 

PRINT ’Longitude limits in files are long( 1 ) ,long(n) 

INPUT "Enter new longitude limits ’; longO.longl 
IF longO<long( 1 ) OR long1>long(n) THEN Inputparameters 
INPUT ’Enter value of Cl ’; cl 
INPUT ’Enter value of C2 ’; c2 
10 INPUT ’Enter value of R ’; R : IF R<0 OR R>1 THEN 10 

20 INPUT ‘Enter value of lambda ’; lb : IF lb<=0 THEN 20 

Functions: 

DEF FNL(tau,Q)=124078!*tau~1.5/Q 

OEF FNQ(Qo,Qs , L)=Qo-R*Qs*(2*lb/L~( 1+2*lb/L)*EXP(-L/lb) ) 

DEF F NFi(taux,Q) = 144. 6*taux*COS( lat/57. 2958 )/Q 

Calculate: 

sum=0 

FOR 1=1 TO n 

IF long(i)<long0 OR long(i)>long1 THEN 30 
IF Qo(i)<=0 THEN h(i)=0 : c$=’*’ : GOTO display 
L=FNL(tau(i),Qo(i) ) 

Q=FNQ(Qo(i) ,Qs(i) , L):IF Q<=0 THEN h(i)=0:c$=’*’:GOTO display 
L=FNL(tau(i),Q) 

Fi=FNFi(taux(i) ,Q) : IF Fi<-.5 THEN Fi=-.5 
h(i)=c1*L/( 1+c2*Fi) 
sum=sum+(depth(i)-h(i) )"2 
c$=’ ’ 
display: 

PRINT long (i );c$ , h ( i ) 

30 NEXT i 

PRINT "Sum : ’; sum 

Choose: 

INPUT ’ Save results( 1 ) ,New f lies ( 2 ) , new parameters(3), Quit ( 4 ) o 

ON o GOTO Saveresults, Inputdata, Inputparameters, quit 

Saveresults: 

INPUT ’Enter name of output file ’; fo$ 

OPEN ■out\*5,fo$ 

FOR i=1 TO n 

WRITE#5,long(i) ,h(i) 

NEXT i 
CLOSE 

GOTO Choose 

Quit: 

END 
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’ Program Spline - fits a cubic spline over a set of (x,y) points. 

' Joaquim Filipe Gaspar - January 1985 - Naval Postgraduate School 



DIM x(200),y(200) , a(200,4) , s(200) : DEFINT i,j,n,e 
DIM a0(200) , a 1(200), a2(200) , a3(200) 

Enterdata: 

INPUT ‘Read disk file(d) or enter data from keyboard(k) *; a$ 
IF a$=’k‘ THEN key ELSE IF a$=‘d‘ THEN disk ELSE Enterdata 

disk: 

INPUT ‘Enter name of file to fit *; f$ 

OPEN *in‘,#1,f$ 
n=1 

WHILE NOT EOF(1) 

INPUT# 1,x(n),y(n) 
n=n+1 

WEND 

CLOSE# 1 

n=n-1 : GOTO Endc 
key: 

INPUT ‘Enter number of points *; n 

FOR i=1 TO n : PRINT i, : INPUT x(i),y(i) : NEXT i 

Endc: 

INPUT ‘Enter end condition ( 1 ,2 or 3) *; endcond 



nm2=n-2 : nm 1=n-1 

dx1=x(2)-x( 1 ) : dy1=(y(2)-y(1))/dx1*6 

defmatrix: 

FOR i=1 TO nm2 
dx2=x(i+2)-x(i+1 ) 
dy2=(y(i+2)-y(i+1 ))/dx2*6 
a(1, 1 )=dx1 : a(1,2)=2*(dx1+dx2) 
a(i,3)=dx2 : a(i,4)=dy2-dy1 
dx1=dx2 : dy1=dy2 
NEXT i 

ON endcond GOTO solve, end2, end3 
end2: 

a( 1 ,2)=a( I ,2)+x(2)-x( I ) 
a(nm2,2)=a(nm2,2)+x(n)-x(nm1 ) 
GOTO solve 
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end3: 

dx1=x(2)-x( 1 ) : dx2=x(3)-x(2) 
a ( 1 ,2)=(dx1+dx2)*(dx1+2*dx2)/dx2 
a( 1 , 3)=(dx2*dx2-dxl*dx1 )/dx2 
dxn2=x(nm1 )-x(nm2) : dxnl=x(n)-x(nm1 ) 
a(nm2, 1 )=(dxn2*dxn2-dxn1 *dxn1 )/dxn2 
a ( nm2 , 2 ) = ( dxn 1 +dxn2 ) * ( dxn 1 +2* dxn2 ) /dxn2 

solve: 

FOR i=2 TO nm2 

a(i,2)=a(i,2)-a(i, 1 )/a(i-1 ,2)*a(i-l ,3) 
a(i,4)=a(i,4)-a(i, 1 )/a(i-1 ,2)*a(i-1 ,4) 

NEXT i 

a(nm2,4)=a(nm2,4)/a(nm2,2) 

FOR i=nm2~1 TO 1 STEP -1 

a(i,4)=(a(i,4)-a(i,3)*a(i+1 ,4))/a(i,2) 

NEXT i 

FOR 1-1 TO nm2 
s(i+1 )=a(i,4) 

NEXT i 

ON endcond GOTO endl 1 ,end22,end33 
endl !: 

s( 1 )=0 : s(n)=0: GOTO coeff 
end22: 

s( 1 )=s(2) : s(n)=s(nml ) : GOTO coeff 
end33: 

s( 1 )=( (dx1+dx2)*s(2)-dx1*s(3) )/dx2 
s(n)=( (dxn2+dxn1 )*s(nm1 )-dxn1*s(nm2))/dxn2 

coeff: 

FOR i=1 TO nml 
dxl=x(i+l )— x( 1 ) 
a3(i)=(s(i+1 )-s(i))/(6*dxl ) 
a2(i)=s(i)/2 

aUi)=(y(i + 1 )-y(i))/dxl-(2*dx1*s(i)+dxl*s(i+1 ))/6 
aO(i)=y(i) 

NEXT i 

calculate: 

INPUT "Enter interval betwen points dx 
PRINT "x-limits on file were : ",x( 1 ),x(n) 

INPUT ‘Enter new x-limits xl ,x2 
INPUT ‘Enter name of output file f$ 

OPEN ‘out‘,#2,f$ 
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FOR i=l TO nml 
xx=x1-x(i) 

WHILE xl<x(i+ 1 ) 

yy=aO( i Ha I ( i )*xx+a2(i )*xx*xx-»-a3( i )*xx*xx«xx 
WRITER, xl,yy 
xx=xx+dx : xl =x1+dx 

WEND 
NEXT i 

xx=x!-x(nml ) 

WHILE xl<=x2 

yy=a0(nm1 Hal (nml Hxx+a2(i)*xx*xx+a3(i)*xx*xx*xx 

WRITE#2,xt,yy 

xx=xx+dx : x1=x1+dx 

WEND 
CLOSE#2 
GOTO Enterdata 
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' Program Graph - plots an Internally defined function or a file from 

disk. 

' Joaquim Filipe Gaspar - Jan 1 985 - Naval Postgraduate School. 

DEFINT a-z : DEFSNG x,y,p : DIM xf(700) ,y(700) ,yf(700) 

pi=3. 141593 :YNF=3E+38 

menudef: 

MENU 1, 0,1 , 'Function' 

MENU 1,1,1 , 'Enter new function' 

MENU 1 ,2,2, 'Use internal function' 

MENU 1,3,1 , 'Input new file from disk' 

MENU 1 ,4,0, 'Use present file' 

MENU 2,0,1 , 'Resolution' 

MENU 2, 1,1, 'Very High' 

MENU 2,2,1, 'High' 

MENU 2,3,2, 'Medium" 

MENU 2, 4,1, -low- 
menu 3,0, 1, 'Scaling' 

MENU 3, 1,2, 'Automatic' 

MENU 3,2, 1 , 'Enter scaling factors* 

MENU 4,0, 1 ,'X-Limits' 

MENU 4, 1,1, 'Enter x-limits* 

MENU 5, 0,1, ‘Grid’ 

MENU 5, 1,2, 'Yes' 

MENU 5,2, 1 , 'No' 

MENU 6,0,1, 'Title* 

MENU 6, 1,1, 'Enter title' 

MENU 7,0, 1, 'Run' 

MENU 7,1 ,0, 'Start graph' 

MENU 7,2,0, 'New graph' 

MENU 7,3,1, 'Stop-List' 

MENU 7, 4,1, 'Quit* 

MENU ON 

res=4 : autoscale=1 : grid=1 

Menuloop: 

IF xrange THEN MENU 7, 1,1 
ON MENU GOSUB Menucheck 
GOTO Menuloop 

Menucheck: 

ON MENU(0)GOSUB Function, Resolution, Scaling, Limits, Grid, Title, Start 

RETURN 
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Function: 

ON MENUO ) GOTO Newfunction,01dfunction,Newfile,01dfile 

Newfunction: 

MENU RESET : STOP 

01 df unction: 
foption=0 

MENU 1,2,2 : IF fflag THEN MENU 1,4, 1 
RETURN 
Newfile: 

INPUT "Enter name of file f$ 

OPEN "in",#1,f$ 

i=0 : yfmin=YNF : yfmax=-YNF 

WHILE NOT EOF(l) 

i=i+1 : INPUT# 1 ,xf(i),yf(i) 

IF yf(i)>yfmax THEN yfmax=yf(i) : imax=i 
IF yf(i)<yfmin THEN yfmin=yf(1) : imin=i 
WEND 
CLOSE# 1 

xfmin=xf(l) : xfmax=xf(i) : nx=i 
f f 1 ag= 1 : foption=1 : MENU 1,2,1: MENU 1 ,4,2 
RETURN 
Oldfile: 

foption=l : MENU 1,2,1 : MENU 1,4,2 
RETURN 

Resolution: 

res=2~(MENU(l)-1) 

FOR i=1 TO 4 : MENU 2,i, 1 : NEXT i 
MENU 2,MENU( 1 ),2 : RETURN 

Scaling: 

IF MENU( 1 )=1 THEN autoscale=1 : MENU 3, 1,2 : MENU 3,2, 1 : RETURN 
autoscale=0 : MENU 3,1,1: MENU 3,2,2 
INPUT "Enter x-axis step "; xstep 
INPUT "Enter minimum y "; ymin 
INPUT "Enter maximum y "; ymax 
INPUT "Enter y-axis step "; ystep 
IF foption THEN yrange=ymax-ymin 
RETURN 

Limits: 

ON foption+1 GOTO limitl ,limit2 
limit I: 

INPUT "Enter minimum x ";xmin : INPUT "Enter maximum x ",* xmax 
xrange=xmax-xmin 

IF xrange=0 THEN limitl ELSE RETURN 
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limitZ: 

PRINT "In present file, minim and max x‘s are : “; xfmin.xfmax 
INPUT “Enter minimum x xmin : IF xmin<xfmin THEN limit2 
INPUT “Enter maximum x “; xmax : IF xmax>xfmax THEN 11mtt2 
xrange=xmax-xmin 

RETURN 



Grid: 

IF MENU( 1 )=1 THEN grid=l : MENU 5, 1 ,2 : MENU 5,2, 1 '.RETURN 

grid=0 : MENU 5,2,2 : MENU 5,1,1 : RETURN 



Title: 

t$=““ : INPUT “Enter title “;t$ : RETURN 



Start: 

IF MENUC l ) = 1 THEN Graph : RETURN 
IF MENU( 1 )=3 THEN MENU RESET : STOP 
IF MENU( 1 )=4 THEN SAVE “Graph" : SYSTEM 
FOR i=1 TO 6 : MENU i,0, 1 : NEXT i : MENU 7,2,0 
CLS : CALL TEXTSIZE(12) : RETURN 



Graph: 

FOR i=1 TO 6 : MENU i , 0 , 0 : NEXT i 
CLS 

ONfoption+1 GOTO graphl ,graph2 

graphl: 

i=0 

IF autosca1e=1 THEN ymax=-YNF : ymin=YNF 
FOR graphx=4 TO 488 STEP res 
i=i+1 

x=xmin+xrange*(graphx-4)/484 

GOSUB Fx 

IF autoscale=0 THEN 30 

IF y(i)>ymax AND ABS(yO)oYNF) THEN ymax=y(i) 
IF y(i)<ymin AND ABS(y(i)<>YNF) THEN ymin=y(i) 
30 NEXT graphx 

yrange=ymax-ymin : ni=i-1 
GOTO setitle 

graph2: 

i=1 

WHILE xf(i)<xmin 
i=i+1 

WEND 

ii=i : xmin=xf(i) 
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i=nx 

WHILE xf(i)>xmax 
j=i- 1 

WEND 

if n=i : xmax=xf(i) : xrange=xmax-xmin 

IF autoscale THEN ymax=yfmax : ymin=yfm1n : yrange=ymax-ymin 

setitle: 

CALL MOVETO( 100,20) : PRINT t$; 

xlabel: 

CALL TEXTSIZEC9) 

LINE(4,4)-(4,284) : LINE(4,284)-(488,284) 

IF autoscale THEN xstep=xrange/4 
FOR x=xmin+xstep TO xmax STEP xstep 
graphx=-484*(xmin-x)/xrange+4 
PSET(graphx,283) 

IF grid=0 THEN 40 
FOR graphy=280 TO 4 STEP -4 
PSETCgraphx, graphy) 

NEXT graphy 

40 CALL MOVETO(graphx-12, 295) 

IF xOxmaxTHEN PRINT x; 

NEXT x 

ylabel: 

IF autoscale THEN ystep=yrange/4 

FOR y=ymin TO ymax STEP ystep 

graphy=280*(ymin-y)/yrange+284 

PSET(5,graphy) 

IF grid=0 THEN 50 
FOR graphx=9 TO 486 STEP 4 
PSET(graphx, graphy) 

NEXT graphx 

50 CALL MOVETO(6,graphy-2) 

PRINT y; 

NEXT y 

Plot: 

ON foption+1 GOTO plotl,plot2 
plotl: 

FOR i=l TO ni : IF (y(i)>ymax OR y(i)<ymin) THEN NEXT i 
gy=280*(ymin-y(1) )/yrange+284 
gx=4+-( i— t )*res : i=0 
FOR graphx=4 TO 488 STEP res 
i=i+1 
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IF (y(i)>ymax OR y(i)<ymin) THEN flagy=0 : GOTO 5 
f1agy=f1agy+1 : graphy=280*(ymin-y(i ) )/yrange+284 
IF f1agy>1 THEN LINE(gx,gy)-(graphx,graphy) 
gx=graphx : gy=graphy 
5 NEXT graphx 
MENU 7,2,1 
RETURN 

plot2: 

gx=-484*(xmin-xf(ii) )/xrange+4 
gy=280*(ymin-yf(ii) )/yrange+284 
FOR i=ii TO if n 

graphx=-484*(xmin-xf (i) )/xrange+4 
graphy=280*(ymin-yf(i) )/yrange+284 
LINE(gx,gy)-(graphx,graphy) 
gx=graphx : gy=graphy 

NEXT l 
MENU 7,2,1 

RETURN 



Fx: 

y(i)=EXP(-x*x) 

RETURN 
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